################################################################################################
#### 运行配置文件
source /public/home/xxf2019/20231121_singleMuti/config/config.sh


################################################################################################
#### RG那一列替换为细胞类型所在的CB那一列
for bam_file in `ls ${input_dir} | grep sorted.bam | grep merge`
do
new_name=`echo ${bam_file} | awk -F"." '{print $1}'`
samtools view -h ${input_dir}/${bam_file} | \
awk 'BEGIN{OFS=FS="\t"} {for (i=1; i<=NF; i++) {if ($i ~ /^CB:Z:/) {gsub(/RG:Z:[^\t]+/, "RG:Z:"substr($i, 6)); break}} print}' | \
samtools view -bS - -o ${input_dir}/${new_name}.reviewRG.bam
samtools index ${input_dir}/${new_name}.reviewRG.bam
done


################################################################################################
#### 基于单细胞的ATAC数据，需要提供peak和bam，每个细胞类型分别鉴定peak对应的motif以及motif在每个细胞的活性程度z-score
#### 参考以下内容写代码
## https://greenleaflab.github.io/chromVAR/articles/Articles/Applications.html#differential-accessibility-and-variability
## 文章：https://www.nature.com/articles/nmeth.4401
genome_type="hg38"
cpu=20

for line in ` cat ${config_path}/cluster_celltype.csv | tr ' ' '_' | tr '&' '-' `
do
clus=`echo ${line} | awk -F, '{print $1}'`
celltype_name=`echo ${line} | awk -F, '{print $2}'`

## 内存需要较大容易溢出，需要检查是否需要重跑
if [ -f ${log_path}/${clus}_${celltype_name}.log ]
then
error_memory=`cat ${log_path}/${clus}_${celltype_name}.log | grep memory | wc -l`
if [ ${error_memory} -gt 0 ]
then 
echo ${line}
echo ${error_memory}
echo " sh ${scripts_dir}/run_chromvar.sh ${clus} ${celltype_name} ${genome_type} ${cpu} ${config_path}" | \
qsub -N ${clus}_"chromvar" -l nodes=1:ppn=${cpu},mem=200gb,walltime=240:00:00 -q student -d ${qsub_log_path}
fi
else
## 之前没用跑过的重新跑
echo ${line}
echo " sh ${scripts_dir}/run_chromvar.sh ${clus} ${celltype_name} ${genome_type} ${cpu} ${config_path}" | \
qsub -N ${clus}_"chromvar" -l nodes=1:ppn=${cpu},mem=200gb,walltime=240:00:00 -q student -d ${qsub_log_path}
fi
done


#### 所有细胞的合并算一个总的
mkdir -p ${output_dir}/combine/
bed_file=`ls ${input_dir} | grep bed | grep cluster | awk '{print input_dir"/"$0}' input_dir=${input_dir} | tr '\n' ' ' `
cat ${bed_file} > ${output_dir}/combine/combine_output.bed
rm -rf ${output_dir}/combine/merged_output.bed
${bedtools} sort -i ${output_dir}/combine/combine_output.bed | \
${bedtools} merge -i - > ${output_dir}/combine/merged_output.bed
rm -rf ${output_dir}/combine/combine_output.bed

## 运行
cpu=10
echo " sh ${scripts_dir}/run_chromvar.combine.sh ${genome_type} ${cpu} ${config_path}" | \
qsub -N "combine_chromvar" -l nodes=1:ppn=40,mem=1000gb,walltime=240:00:00 -q smp -d ${qsub_log_path}


################################################################################################
#### 1、计算TF的活性和基因表达的相关性
#### 2、根据peak和基因的关系以及peak和TF的关系，计算TF在该基因的激活程度占TF总激活程度的比例
cor_gene_peak=0.5
cpu=5

for line in ` cat ${config_path}/cluster_celltype.csv | tr ' ' '_' | tr '&' '-' `
do
cluster=`echo ${line} | awk -F, '{print $1}'`
celltype_name=`echo ${line} | awk -F, '{print $2}'`

echo ${line}
echo " sh ${scripts_dir}/tf_gene.sh ${cluster} ${cor_gene_peak} ${celltype_name} ${cpu} ${config_path}" | \
qsub -N ${cluster}_"tf_gene" -l nodes=1:ppn=16,mem=100gb,walltime=240:00:00 -q smp -d ${qsub_log_path}

done
